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Abstract 


A new formulation is presented for a variational calculation of A^-body sys¬ 
tems on a correlated Gaussian basis with arbitrary angular momenta. The 
rotational motion of the system is described with a single spherical harmonic 
of the total angular momentum L, and thereby needs no explicit coupling of 
partial waves between particles. A simple generating function for the corre¬ 
lated Gaussian is exploited to derive the matrix elements. The formulation 
is applied to various Coulomb three-body systems such as e“e“e''', ttfu, tdfi, 
and ae~e~ up to L = 4 in order to show its usefulness and versatility. A 
stochastic selection of the basis functions yields good results for various an¬ 
gular momentum states. 
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I. INTRODUCTION 


We have recently shown that the stochastic selection of the correlated Ganssian 
basis fnnctions leads to precise variational solntions for diverse fermionic and bosonic nona- 
diabatic N = 2 — 7-body systems. The trial fnnction in a variational approach must be 
flexible enough to describe the full variety of correlations between the particles. The corre¬ 
lation is conveniently represented by a correlation factor, F = Ylf^j fij, and thus the trial 
function is often chosen to be of this form. The trial function of Hylleraas type, which is 
used in atomic and molecular physics, approximates as a linear combination of expo¬ 


nentials, exp[—ttj 
exp [—a 
n*v,exp[-a 


. r, — r. 


Ti — Tjl]. If fij is approximated as a linear combination of Gaussians, 
the iV-particle trial function then contains a product of these Gaussians: 




Ti — r. 


= exp[- aij 


(U - rjf] 


The above form of the correlated function describes the motion with the orbital angular 
momentum L = 0 only. It is obviously important to extend the correlated function to the case 
with L > 0, and in fact there is an increasing interest in hnding a precise solution with L >2 
for the Goulomb three-body systems [0 . The standard way to describe the rotational motion 
is to vectorially couple the solid spherical harmonics of the relative coordinates |]^,^. Each of 
the solid spherical harmonics carries the partial wave for the corresponding relative motion. 
Since these partial waves are in general not good quantum numbers, several sets of partial 
waves are in general necessary for a realistic description of the motion, particularly in a 
nuclear system [§,|^ . Another way to incorporate the angular dependence is to use Gartesian 
Gaussian functions. The calculation of matrix elements becomes complicated in both cases, 
especially when the number of particles increases or high angular momenta are involved. 
A quite different way to introduce the angular part, proposed in [Q, involves no partial 
wave decomposition for each relative motion but attempts to determine a vector which has 
closest relevance to the rotational motion. The vector is defined as a linear combination of 
the relative coordinates and their coefficients can be treated as variational parameters. The 
formulation using this new angular part is entirely free from the complexity involved in the 
angular momentum coupling. 

The purpose of this paper is to exploit this formulation in greater detail, to present 
the formulas needed for A^-body system interacting via central force, and to test its utility 
by applying to Goulomb three-body problems such as the muonic molecule and the helium 
atom. Several authors have investigated both the ground and excited states of these systems, 
for example in |^,]T-]G] for the muonic molecule and in 


15,161 for the helium atom. A 


comparison with the solutions known in literature will be useful to judge the utility of 
the present formulation. Although any spherically symmetric orbital functions can be used 
together with the angular part introduced here, the correlated Gaussian has the advantage 
of simplicity in the coordinate transformation. 

In section II we introduce the correlated Gaussian with the angular momentum de¬ 
pendence and show that the matrix elements can easily be evaluated with the use of the 
generating function of the correlated Gaussian. No problem arises as to the center-of-mass 
motion because the formalism does not include any dependence on the center-of-mass vari¬ 
able. In section III some of the results for the Goulomb three-body system are presented 
within the new formulation for the angular dependence and compared to those available in 
literature. A summary is given in section IV. 
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II. FORMALISM 


A. The correlated Gaussians 

Any square-integrable function with angular momentum Im can be approximated, to any 
desired accuracy, by a linear combination of nodeless harmonic-oscillator functions (Gaus¬ 
sians) of continuous size parameter a 

rim(r) ~ with = |r|'y';^(f). (1) 

A generalization of this to A-particle systems contains a product of the Gaussians 
exp[—J2f<jOiij{ri — The product can be conveniently expressed in terms of a set of 

(iV — 1) independent relative coordinates x, (xi, ...,xjv-i), instead of N{N — l)/2 interpar¬ 
ticle distance vectors (rj — Fj). By a set of relative coordinates we mean the one in which 
the intrinsic kinetic energy operator takes the form with reduced masses /ij. 

Even with this condition there are a number of possible sets of relative coordinates. One can 
choose, however, any one of the sets as x because each set of relative coordinates is obtained 
from any other set of relative coordinates by an appropriate [N — 1) x {N — 1) matrix 
T. The relative coordinates are assumed to be normalized in such a way that the volume 
element remains unchanged under the coordinate transformation, which requires that the 
determinant of T is unity. An iV-particle basis function, the so-called correlated Gaussian, 
then looks like 


171 


= e (2) 

where A is an {N — 1) x [N — 1) positive-definite, symmetric matrix containing N{N — l)/2 
nonlinear parameters, specihc to each basis element, and the quadratic form, xAx, involves 
scalar products of the Gartesian vectors: xAx = AjXj ■ x^. 

The function in Eq- (2), which represents the angular part of the wave function 

with the total orbital angular momentum L and its projection M, is a generalization of y 
in Eq. (1). Usually it is chosen as a vector-coupled product of solid spherical harmonics of 
the relative coordinates 

^Lm(x) = [[[yh{^l)yh{^2)W2yi3i.^3)]Li23--]LM ’ (3) 

where the square bracket stands for the coupling of angular momenta. Each relative motion 
has a definite angular momentum in Eq. (3). Since the set of angular momenta itself is 
not a conserved quantity, it may be important to include several sets of angular momenta 
(/l, I 2 , ..., In-U Li2, Li23, ...) for a realistic description. This is the case especially in nuclear 
few-body problems P, pT|JT8[| . It is also noted that a faster convergence is in general obtained 
by allowing the use of different sets of relative coordinates together with suitable sets of 
angular momenta. From the fact that 6 '£m(x) can be expressed by different partial wave 
decompositions in different relative coordinate systems, one can conclude that the usage of 
partial waves may not be so important after all. Besides, the various possible partial wave 
contributions increase the basis dimension. Moreover, the calculation of matrix elements 
for this choice of 6lm{^) sooner or later becomes too complicated. This choice is therefore 
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apparently inconvenient especially as the number of particles increases and/or the different 
sets of relative coordinates are employed. 

As proposed in |^, this difficulty can be avoided by adopting a different choice for 

N-l 

6 'lm(x) = r]KLM{u,x) = |vP^+^1Am(v), with v = ^ = ux. (4) 

i=l 


Only the total orbital angular momentum, which is (at least approximately) a good quantum 
number in most cases, appears in this expression. The real vector u = (mi, ...,MAr_i) dehnes 
a global vector, v, a linear combination of the relative coordinates, and the wave function of 
the system is expanded in terms of its angle v. The vector u may be considered a variational 
parameter and one may try to minimize the energy functional with respect to it. The energy 
minimization then amounts to hnding the most suitable angle or a linear combination of 
angles. The continuity of the parameter u can be more advantageous in a variational calcula¬ 
tion than the discrete nature of the set of the angular momenta (/i, h ,..., In-i, L 12 , T 123 ,...) 
because the change of the energy functional can be continuously seen in the former case. 
The factor of plays an important role in improving the short-range behavior of the 

A remarkable advantage of this form of 


IS 


wave function, e.g., the Coulomb cusp ratio 

is that the calculation of matrix elements becomes much simpler than in the former 
case because the coupling of (A^ —1) angular momenta is completely avoided. 

The two forms of 6 'lm(x) are in fact closely related to each other. It is easy to see that any 
of the functions of Eq. (4) is a linear combination of the terms of Eq. (3), each multiplied by 
an appropriate monomial of the variables, xi^,..., x^r_i^. For example, it takes a particularly 
simple form for the three-body system, i.e., for the vector v = UiXi -|- U2X2 


,|2K+Ly^ 


lmIv, 


= 


E 


2fci+/i 2/C2+^2 

Ui 112 


fcj li ,/c2 12 —^ 


X 


(2A:)!!(2A: + 2L + 1)!! 


( 2 li 1) (2/2 + 1) 


{2ki)\\{2h + 2h + l)!!(2fc2)!!(2A;2 + 2/2 + 1)!' 


2L + 1 


(/ 1 O/ 2 OILO) 


X xd 


X2 


yh{^i)y 12(^2) 


J LM 


(5) 


Table I lists possible sets of fci, /i, ^ 2 , and I 2 values for small K and L values. In the case of 
K = 0 both ki and ^2 are limited to zero and only the stretched coupling, namely / 1 -I -/2 = L, 
is allowed. With an increasing K value the possible values of partial waves li and I 2 increase 
including the case of non-stretched coupling. To increase K is thus one way to include 
higher partial waves in the calculation. Note, however, that even with K = 0 additional 
and important partial wave contribution comes from the cross term of the exponential part 
of the correlated Gaussian if A is not diagonal. Conversely any two-variable functions of 
Eq. (3) with natural parity, i.e., (—l)^i +^2 _ pg expressed in terms of a linear 

combination of the terms, |vp'^+^y/M(v), by using some appropriate sets of u values, each 
multiplied by a monomial of degree I 1 + I 2 — 2K — L in xi^ and X 2 ^. See Appendix for the 
details. A generalization of the argument in Appendix will lead to a conclusion that some 
functions of Eq. (3) may be given in terms of a linear combination of Eq. (4), each multiplied 
by the terms such as and (xj ■ xj). Namely, the rotational property of the function of Eq. 
(3) may be represented by combinations of simple forms of Eq. (4). Therefore, if one can 
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calculate the matrix elements using defined in Eq. (4), then those with the previous 

form of of certain class can be obtained readily. 

The correlated Gaussian we proposed is thus given by 

fxLMiu, A x) = t]klm{u, x)e“5^^^. (6) 

A useful property of the function / is its form-invariance with respect to the transformation 
of the coordinates x to any other set of coordinates y, that is, for x = Ty, 

a, x) = fxLMiu , A', y), (7) 

with 


A' = TAT, u = fu. (8) 

Here T is the transpose of T. As will be seen later, this property is fully exploited in evalu¬ 
ating the matrix elements. It is necessary to impose a proper permutation symmetry on the 
basis states for the system of identical particles. The symmetry requirement causes a linear 
transformation of the coordinates, and it can easily be incorporated in the present formu¬ 
lation thanks to the form-invariance of the correlated Gaussian mentioned above. It would 
be rather complicated to construct the symmetry-adapted basis states using the angular 
function of Eq. (3). 

The calculation of the matrix elements becomes simpler if one uses a generating function 
of the correlated Gaussian. In fact, the following function g is found to be convenient to 
generate the function / |I|: 

IkimA, a, x) = —— J dtYiMit) ( A, x) j , (9) 

where 

g{a, t; u, A, x) = 


_ 4:7i{2n + l)\ 

“ 2^n!(2n + 2/ + l)!!‘ 

Here t is a unit vector. Equation (9) is easily proved by using the simple formula 
(a-b)" = la^bl" ^ Br,i ^ YiUAYimAT- 


n,l>0 m= — l 

2n+/ = fc 


(11) 


( 12 ) 


The correlated Gaussian basis of Eq. (6) has parity (—1)^. To construct a function with 
parity (—the angular part of Eq. (6) must be slightly modified, e.g., to 

OlmA) = [^7KL(^^,x)?7ol(^^^x)]iM• (13) 

In this case the generating function g of Eq. (10) must be modified to include another factor 
a't' ■ (u'x.) . Equation (9) is then extended to 


[r]KL{u,AVoi{u',A]LMe 

^ ^ ^dtdt'flA(t)Ti(t') 


BklBqi 


/ ^ _g—ixylx-|-at'(?ix)-|-a't'(ii'x) 




Oi=ot'=0 

t=t' = l 


(14) 
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B. Calculation of the matrix elements 


In this subsection we will give the details of the method of calculating the matrix elements 
between the basis function of Eq. (6). One can already see some useful formulas for L = 0 
motion in The aim here is to demonstrate that the matrix elements for arbitrary L 
can be obtained as simply as those for L = 0. As will be seen later, one generally needs 
the matrix elements of the two functions which are expressed in terms of different sets of 
relative coordinates. 

The overlap matrix element is obtained as 


{fKLM{u,A,y)\fK'LM{v,B,^)) = (/ 

-D- [ [ 

okl^jk'l J J 

( (27r) 


ijiVl / \J Y 


N-l \ 3/2 


da'^da"^' detO 
with an abbreviation 


+p'ci'^ +qaa'tt' 


(16) 


a = a' =0 
t=t' = l 


K = 2K + L, k' = 2K' + L, 
where the property of Eq. (7) is used and where 

C = A + B', p =-uC~^u, p' = ^v'C~^v', q = ^(uC~^v'+ v'C~^v\ . 


(16) 

(17) 


Here we used the familiar formula of the 3n-dimensional Gaussian integration 

where A is an n x n symmetric matrix and T = (Ti, ■ • T„) is a row vector comprising three- 
dimensional vectors Tj. Differentiating with respect to a and a', followed by a = a' = 0, 
and integrating over the angles of t and t' leads us to 


{fKLM{u, A, y)lfK'LM(v, B, x)) 

\N-1\ 3/2 mm{K,K') 

5o iK-ny.iK' -n)\{L + 2n)\ 


1 

BklBk'l 


■{2^ 

detG 


pK-np/K'-n^L+2n 


B. 


nL' 


(19) 


The matrix element of the kinetic energy operator is calculated in a similar way. Ex¬ 
pressing the intrinsic kinetic energy operator in terms of the set of coordinates y of the bra 
side, we obtain 

TV—1 2 A^—1 2 

{fKLM{u,A,y)\ ^ ^\fK'LM{v, B,x)) = {f kLm{u, A,y)\ ^ ^\fK'LM{v',B',y)) 

i=l i=l 


1 


— J y'dtdt'lLM(t)*yLM(t')f {g{a,t,u,A,y)\ 


Bkt.B 


KL-Dr'L 


X 


STrAH' — a'^v'Av' + 2a't' ■ (v'AB'y) — yB'AB'y \g{a', t', v', B', y)) 


ctc = a'=0 
t=t' = l 


. ( 20 ) 
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Here A is an [N — 1) x [N — 1) diagonal matrix 


0 - 0 \ 


A = 


2/12 


V 


0 




( 21 ) 


2/ijv- 


-J 


The integration over y in Eq. (20) can be done by using the formula of Eq. (18) as 
follows: 


dx 

3/2 

2 ' 
dx 

3/2 


'A. 

da 


dx e-h^MT+aV)> 


Q :=0 




2 A /rfx e-4 

da J 


5x(yl+aS)x+Tx \ 

/ 0=0 


) (3TrA-^H + TA-^HA-'T)e 


-1 D /l-l^r^„iTA-lT 


Using Eqs. (22) and (23) in Eq. (20) leads us to 


V-l 2 

(/klm(m, A, y)| ■^\fK'LM{v,B,x.)) 

i=l 


-JJ didi'YLM{iyYLM{i')\ 


BklBk'l 

\N-l \ 3/2 


d«+«' 

da'^da'^' 


X 


(2ir) 


detC 


R + Pa^ + P'a^ + Qaa't ■ t' 


^po^+p'o'^+goo't't' 


q:=q:^ = 0 
t=t' = l 


where 


R = 3TtAB'C-^A, P = -uC-^B'AB'C-^u, P' = -UC'-^AAAC'-^U, 
Q = uC-^B'AAC-\' + PC-^AAB'C-\. 

The matrix element of the kinetic energy operator is hnally obtained as 

N-l 2 

(/xlm(m, A, y)| Y : 5 ^l/x'LM('y,-B,x)) 

i=l 


’(2t 


Bp't^Brc't, \ dctCx 


TV—1\3/2 mm{K,K') 

K,\k,'\ Y, 


n=0 


X 


Ppp'q + {K — n)Pp'q + {K' — n)pP'q + (L + 2n)pp'Q 


pK—n—1 pf qL-\-2n—l 

^(K- n)\(K' - n)!(L + 2n)! 


( 22 ) 


(23) 


(24) 


(25) 


(26) 
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The present formulation does not have any problem arising from the center-of-mass motion 
as discussed in 

We will show the evaluation of the potential energy matrix element by assuming that the 
potential Vij is a function of the distance |rj —r^l only. There is at least one set of coordinates, 
say, z, in which (r^ — rj) can be chosen to be zi. It is then convenient to calculate the matrix 
element of the potential by transforming both of the correlated Gaussians on bra and ket 
sides to those expressed in this set of coordinates: 


{fKLM{u, A, y)\Vij\fK'LM{v, B, x)) — {fxLMiu', A\ 7,)\V {zi)\f k'LmW, B\ z)) 

—- - - f f didi'YLM{i)*YLM{i') 

^KLJ^K’L J J 


da'^da"^' 


{g{a, t, u\ W, z)\V{zi)\g{a', t', n', B\ z)) 


a. = a.' =0 
t = t' = l 


(27) 


where A' and u' are dehned by Eq. ( 8 ) with the matrix T corresponding to the y ^ z 
transformation y = Tz, while B' and v' are obtained similarly by the matrix corresponding 
to the X —z transformation. 

To perform the integration over z in Eq. (27) we introduce the short-hand notation as 
follows: 


u' = {u[,U 2 ,-- ■,Un-i) = v' = = Kw), 

/ C 7i ... 7Ar_2 \ 

7i 


A' + B' = 


V 7V-2 


(28) 


/ 


By introducing 7 = ( 71 , • • •, 'yN- 2 ) and z = (zi, • • •, z^^i) = (zi, w), the z integration of Eq. 
(27) becomes 


{g{a,t,u',A',z)\V{zi)\g{a,t',v',B',z)) = JJ dzidw V{zi) 

xexp^- czi'^ + ■ zi-wTw -|- at ■ (cDw) -|- cit! ■ (yw) — zi ■ (qw) ^ 


2 

\N-2 \ 3/2 


j ^va'^+pW^+lvaa't-t’ j (^ 29 ) 

79„'|Ar-2\ 3/2 

^ + AVt>i). 

detr J Jo 


where 


Pv = ^o;, = ^yr ^y, 

X = u[-^ fyT-^o; + (hT’S j , 


(lv = \ ^y + yT . 

A'-wl-ihr-y + xr-x), 
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( 30 ) 










and iQ^x) = smhx/x. Substituting Eq. (29) into Eq. (27), followed by a power series expan¬ 
sion of Iq^x), and doing the needed operation leads us to 


{fKLMiu, A, y)lV(zi)lfK'LM(v, B, x)) 
1 fi27r)^-^V^" 


Bkt.B 


KLJ^K’L 


detr 


K+K'+L m m-i m\n{K-i,K'-j) 

H /(2m + 2,c-7r“S) 

m=0 i=0 j=0 n=ma.x(0, ^-^-^-^ ) 


X 


{2m + l)\i\j\{m — i — jy.{K — i — n)\{K' — j — n)\{L + 2n — m + i + j)\ 


B., 


nLi 


with 


poo 

I{n,a) =4:71 dzi V{zi)zi 
Jo 




(31) 


(32) 


The integral of Eq. (32) becomes elementary for the Coulomb potential. The calculation of 
the mean distance, the root mean square distance or the mean inverse distance between the 
particles can easily be done by putting V{z) = z^ z^ oy Xjz m. the above integral. 

The calculation of the matrix elements described above is much simpler than the case 
where the function is decomposed into partial waves of the relative coordinates as 

in Eq. (3). In fact, in that latter case one has to integrate over the angles of the relative 
coordinates and one has to cope with the angular momentum algebra. We note, however, 
that the calculation of the matrix element of the latter type poses no problem if the function 
OlMl (x) of Eq- (3) is expressed as a linear combination of the terms of Eq. (4) with appropri¬ 
ate M-vectors. It is appealing that the present formalism does not require any modihcation 
with respect to an increasing N. 

All the matrix elements can be given in a closed analytic form and the numerical evalu¬ 
ation of the matrix elements as a function of the nonlinear parameters is therefore straight¬ 
forward. The values, k = 2K+L and k' = 2K'+L, are usually small in practical cases and 
the sum in Eqs. (19), (26), and (31) is limited to just a few terms. 

The partial derivative of matrix elements with respect to variational parameters may 
sometimes be useful when one searches for an optimal set of parameters. Since the de¬ 
pendence on the parameters is explicitly given for the matrix elements of the correlated 
Gaussians, it would not be difficult to derive the expression for the derivative. As stated in 
the previous section, the w-vector is considered a variational parameter which dehnes the 
most suitable global vector to describe the rotational motion. The calculation of the deriva¬ 
tive of the matrix elements with respect to is particularly simple because of the simple 
structure of the w-dependence. It can also be calculated by using an equation analogous to 
Eq. (14) because the derivative of tjklm with respect to Ui is expressed as a tensor product 
of two ? 7 ’s as follows: 


^Vklm{u,x) = {2K + 2L + 1)^^^^^ 


?7ifL_l(M,X)Xj 


LM 


2Kx 


I L + 1 
2 L + 1 L 


r]K-iL+i{u,yi)yii 


\LM 


( 33 ) 


Note that x* is equal to An/?)rjQi{u', x) with u' being dehned by -u) = 5ij. 
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III. SPECIFIC EXAMPLES: COULOMB THREE-BODY SYSTEM 


This section is devoted to show the bound-state solutions of Coulomb three-body prob¬ 
lems by using the method described above. To test the method we consider the positronium 
negative ion e“e“e^, the muonic molecules of ttfi~ and tdfi~, and the helium atom ae~e~. 
These systems are chosen because they cover a wide range of mass ratio of the constituent 
particles and include some rotational bound states: The positronium negative ion consists 
of particles of identical mass, which makes it difficult to treat the system in an adiabatic 
approach, the mass of t or d in the muonic molecule is a few ten times heavier than the 
mass of p, and in the helium atom the mass of the heavy particle is several thousand times 
larger than the e~ mass, which is the realm of the adiabatic approach. One can in principle 
consider other systems with the method, but there are already highly accurate calculations 
available for the above systems. We assume that the orbital wave function is antisymmet¬ 
ric with respect to the interchange of the two identical particles except for the case of the 
spin-singlet states of the helium atom. 

The trial wave function is chosen to be a linear combination of the correlated Gaussian 
of type (6). Without loss of generality the vector u can be set to satisfy uf + = 1. Each 

basis function thus contains at most four nonlinear parameters, three of which come from 
the matrix A. To assure the positive definiteness of the correlated Gaussians, A can in 
general be expressed as A = GDG [^], where G is an (N — l) x (N — l) orthogonal matrix 
containing (A —l)(iV —2)/2 parameters and H is a diagonal matrix, Dij = diSij, including 
(N — l) positive parameters di. Although no restriction on the parameters of the matrix G 
is in principle necessary, it is advisable to avoid too many variables if possible. The most 
naive choice would be to take G as a unit matrix, which is equivalent to using only a single 
set of coordinates x, and then to try to reach convergence by including successively higher 
partial waves implied by an appropriate choice of K and u values. Many examples show 
||I|J8|JT^ , however, that this type of single channel calculations does not work well especially 
in those systems where the adiabatic approximation of the motion is questionable because 
the various type^of correlated motion become important. Instead the matrix A may be 
chosen as A = T~^DT~^ , where T is such a matrix that connects the coordinates x to an 
arbitrary set of relative coordinates y, x = Ty. In this case xAx becomes yDy = diyf. 

The Gaussians with this choice of A together with the angular part of Eq. (3) have recently 
been shown to give a precise solution for few-nucleon systems with realistic nucleon-nucleon 
potentials including noncentral components |10|. Since our interest is to test the utility of 


Eq. ( 6 ), we restricted the choice of A and a maximum value of K to the following three 
simple cases: 

(i) Amax = 0 and A = T~^DT~^, where T represents a special matrix connecting to three 
sets of relative coordinates (the so-called rearrangement channels) and D is diagonal. 
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A m a x = 1 and the choice of A is the same as in case (i) 


111 


A m a x = 0 but A is an arbitrary positive-definite, symmetric matrix. 


In case (hi) the matrix A can be parametrized as A = GDG with an orthogonal matrix G. 
As was already mentioned in the previous section, the angular part with A = 0 describes 
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the stretched conhguration and therefore allows such limited angular correlations in case (i). 
The case (ii) is an extension of this case to include the non-stretched coupling. The angular 
correlation between the particles is taken into account in case (iii) by the cross term of the 
exponential part of the correlated Gaussian. 

Our primary purpose is to investigate if our basis functions with the new angular part can 
universally yield satisfactory results of the same quality independently of the system. Results 
of the calculation presented below are thus limited up to the maximum basis dimension 
/C = 200 for all the systems. With a larger basis size one can dehnitely obtain much more 
accurate results but this is beyond our purpose. Even with /C = 200 dimension there are a 
large number of nonlinear parameters. A complete optimization is in principle superior to 
any other method but it is neither practical nor possible at present. Instead we select the 
basis elements in the trial function stepwise according to the stochastic variational method 
(SVM) ||l| J^J^Jl8|JT0[| . The SVM attempts to select the most appropriate basis elements in 
a trial and error procedure: Various randomly generated candidates for the basis element are 
tested and the usefulness of these states are judged by their contribution to the energy of 
the system. The SVM has proved to provide a precise solution for various few-body systems 
with a reasonable computational effort. The rehning procedure employed in to tune 
the nonlinear parameters was effective to reach the solution of high quality for three- and 
four-nucleon systems interacting with realistic nuclear potentials and was successfully used 
in the present paper as well. The basis selection in the SVM can be done by exploiting the 
special form of the Hamiltonian matrix as shown in Ref. [|I| and it does not carry heavy 
computational loads. There are several sophisticated optimization strategies in quantum 
chemistry IBQ including, for example, the random tempering PJT3[|. The main difference 
between the SVM and other random optimization methods such as the random tempering 
is that the SVM employs the step-by-step procedure to build up a basis set, optimizing the 
nonlinear parameters of the basis states with respect to each other. The advantage of using 
the SVM here is that it seems to be quite suitable to hnd the optimal u vector. 


The results for the positronium negative ion are compared to other calculations [22-27 


in Table II. The energy obtained at the dimension /C is also shown. The trial function of case 
(i) does not produce good result, which is not surprising because the interparticle correlation 
is poorly represented in this case. The result of case (ii) shows a signihcant improvement over 
the case (i), conhrming the importance of the polynomial part of type (4). The energy and the 
root mean square radius in this case are slightly better than our previous calculation |l| . This 
improvement is due to the rehnement in the optimization of the nonlinear parameters, which 
were not employed before. A full correlated Gaussian basis of case (iii), though no polynomial 
part is employed, gives even better results which are comparable to the results obtained in 


the similar basis dimension by using the generalized Hylleraas-type wave functions 22,23 


Gompared to the correlation-function hyperspherical harmonic (GFHH) method |p5 


our 


correlated Gaussian basis seems to give a better solution. The table also lists the calculated 
average distances between the two electrons and between the positron and an electron. They 
are in reasonable agreement with the most precise values obtained by using the correlated 
exponential (GE) functions or the Hylleraas-type functions . 


Table HI lists the results for the lowest L = 0 —3 states of the ttfi molecule. This system is 
analogous to the negative positronium ion though the mass ratio of the constituent particles 
is considerably different. The results with the cases of (ii) and (iii) are shown and the quality 
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of them is similar to what is mentioned in the case of the negative positronium ion. Our 
calculation for the S and P states reproduces the hrst seven digits of Ref. which uses 


the correlated Slater-type geminals or interparticle CE functions. For the D state our result 
is in good agreement with the extensive calculation of |T^ using JC = 2250 basis functions 
similar to [^. Both of them employ the bipolar harmonics of the stretched coupling to 
describe the rotational motion. The calculation of case (iii) uses R"max = 0 so that the basis 
function in this case also employs only the stretched coupling. The fact that we have got 
rather accurate results in a small basis suggests that the usage of Eq. (6), particularly its 
angular part is very useful. We conhrm that the F state is bound, in agreement with the 
other calculation 


111, but have found no bound G state. 


We next show in Table IV the results of calculation for the dtp molecule. The basis 
functions of case (iii) again give better energies than those of case (ii). In fact they reproduce 
the hrst six digits of the most precise variational calculations for L = 0 — 2 states. 

The Gaussian basis similar to case (ii) is employed in the GBCRC calculation of |^, where 
the angular part is, however, represented by the successive coupling of type (3). The fact 
that the D state energy of our calculation with JC = 200 becomes slightly lower than that of 
the GBGRG calculation with JC = 1566 conhrms that a careful optimization of the nonlinear 
parameters is very important and moreover that the angular function of Eq. (4) is really 
useful. In Ref. the optimized energy at smaller JC values is given for the S and P 


states. Our energies are slightly lower than theirs at the same dimension, which indicates 
the effectiveness of the basis optimization of the SVM. 

Finally we present in Table V the energies of the helium atom which are obtained with 
cases (ii) and (iii). We used the hnite mass of the a-particle, rua = 7294.2618241me. The 
basis functions of case (iii) give lower energies than those of (ii) for the S and P states but 
produce considerably higher energies for the states with higher L values. This suggests that 
the explicit introduction of the non-stretched coupling is more favorable to describe in our 
formulation the rotational motion of highly asymmetric states. In the table the energies of 
the S and P states are compared to the result of Ref. 


15 


where the bound states of the 
atomic and mesomolecular systems were studied up to L = 2 by using the GE functions 
of the interparticle distances. The angular part is expressed by the bipolar harmonics. The 
agreement is fair. To include the comparison with the combined conhguration interaction 
Hylleraas method calculation |^, we repeated the calculation for the D and F states as¬ 
suming the inhnite mass for the a-particle. We see that our basis functions reproduce the 
energies of both states up to the hrst six-seven digits. The advantage of our method is that 
no special care is needed to treat states with high angular momenta. 


IV. SUMMARY 

We have presented a formulation for the correlated Gaussian basis which can describe the 
orbital motion with an arbitrary angular momentum. Instead of the well-known successive 
couplings of the spherical harmonics to the total orbital angular momentum the angular 
part of the basis functions in this formulation is uniquely represented by a single spherical 
harmonic of a global vector which is a linear combination of the relative coordinates. We 
have shown that this type of the correlated Gaussians is simply derived from the generating 
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function which is invariant with respect to the transformation of the relative coordinates. 
This property makes it much simpler to calculate the matrix elements and in fact the for¬ 
mulas for the matrix elements have been derived in a compact form for a general system of 
7V-particles interacting with the central forces. They can be applied universally to diverse 
systems independently of the number of particles and of the total angular momentum. The 
extension to include the noncentral forces is straightforward. It should be noted that the 
present formulation for the angular part can be applied to other orbital functions as well as 
the correlated Gaussian. 

The formulation has been applied to obtain the bound states of various Coulomb three- 
body systems which cover a wide range of mass ratio, namely, e“e“e+, ttfi, tdfx, and ae~e~. 
The basis functions have been set up stepwise in a trial and error procedure using the 
stochastic variational method. The coefficients dehning the global vector can be varied con¬ 
tinuously to minimize the energy. The nonlinear parameters of the correlated Gaussians are 
parametrized by three different options. One of the options is to include only the so-called 
rearrangement channels with the angular functions which describe the non-stretched cou¬ 
pling. This model has already given fairly good results. Compared to the restricted case in 
which only the stretched coupling is included, this option presents much better results. The 
other option is to use the full matrices to include the correlation between the particles. This 
model, even without use of the polynomial part, has been better than the former model 
in most cases but the high L states of the helium atom and, with the basis size of 200, 
reproduced the hrst six—seven digits of the total binding energies of almost all the systems 
mentioned above. We have tested the bound states up to L = 4 and found that no difficulty 
arises from the treatment of high angular momentum states in our formulation. We have 
conhrmed that the correlated Gaussian of type (6) works very nicely as a basis function for 
bound-state solution for diverse systems. Although we suggested the method to treat the 
unnatural parity case, a more unihed representation for both natural and unnatural parity 
cases should be developed. 

Finally we summarize some merits of our method in the following. 

(i) No partial wave expansion is needed, and thus no problems in angular momentum 
coupling arise. 

(ii) Universality of the scheme. One needs to introduce no change, for example, when treat¬ 
ing a larger system of Wparticles or describing states with high angular momentum. 
There is no center-of-mass motion problem. 

(hi) Invariance with respect to the coordinate transformation. The form of basis states is 
kept invariant under the transformation. It is easy to construct the symmetry-adapted 
basis states. 

(iv) Fully analytical calculational scheme. The matrix elements are evaluated as simply as 
those for L = 0. The dependence of the matrix elements on the nonlinear parameters 
has simple structure. 
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APPENDIX 

We prove by induction that the vector-coupled function of two solid spherical harmonics 

With (-!)'■+'> = (-l)i (34) 

can be expressed in terms of a linear combination of terms 

|xy>’'|x2p>’"|v|2»+WiM(v). (35) 

where the degree of the function is given by 2pi + 2p2 + 2q + L = li + I 2 , and the vector v 
of each term is given by v = wiXi -|- ^ 2 X 2 with appropriate coefficients Ui and U 2 . 

First we prove that the statement is true for a special case of li + I 2 = L, namely for the 
lowest order terms for a given L. As Eq. (5) shows, |v|^Wm(v) consists of L -|- 1 terms of 
A’/(xi)3^l_z(x2 ) (/ = 0,1, • • •,//), each multiplied by u\u 2 ~^. By using L -|- 1 mutually 

different Ui values, in Eq. (5) and keeping U 2 an arbitrary constant, 

it is possible to pick up a particular term 3 ^i(xi) 3 ^i_;(x 2 ) through a linear combination 

of Eaii Co|v„|^WM(Va) with V„ = Mi“^xi -|- M 2 X 2 , where Cq, should satisfy the following 
equation 

L+1 

^(m^“^)™c« = JmZ, for m = 0,1, •••,/.. (36) 

(y.=l 

The solution of the above linear equation for Cq, is given by Vandermonde’s determinant. 
Thus our assertion is proved. 

Next we assume that the statement holds for all the cases of I 1 +I 2 2{K — 1) + L. We 

note that all the terms of yK+i{^i)yK+L-i{'^ 2 ) (^ = 0,1, ■ ■ ■,//) with degree 2K + L 

appear in the expansion of Eq. (5) together with the respective coefficients The 

rest of the terms in the expansion can, by assumption, be expressed in the form of Eq. (35) 
because the vector-coupled part of each term has smaller degree than 2K + L. Using L -|- 1 
different Ui values similarly to the above case of Eq. (36) enables us to express each of the 
terms, Vx+z(xi)Vx+l-z(x 2 ) , in terms of a linear combination of the functions of Eq. 

L J LKd 

(35). This completes the proof. 
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TABLES 


TABLE I. Possible sets of partial waves contained in a single spherical harmonic with low L and 
K values for a three-body system. See Eq. (5). Listed below are only the cases of 2ki + li < 2/c2 + ^ 2 - 


L 

K 

ki 

h 

k2 

h 

0 

0 

0 

0 

0 

0 


1 

0 

0 

1 

0 



0 

1 

0 

1 


2 

0 

0 

2 

0 



0 

1 

1 

1 



0 

2 

0 

2 



1 

0 

1 

0 

1 

0 

0 

0 

0 

1 


1 

0 

0 

1 

1 



0 

1 

0 

2 



0 

1 

1 

0 


2 

0 

0 

2 

1 



0 

1 

1 

2 



0 

1 

2 

0 



0 

2 

0 

3 



0 

2 

1 

1 



1 

0 

1 

1 

2 

0 

0 

0 

0 

2 



0 

1 

0 

1 


1 

0 

0 

1 

2 



0 

1 

0 

3 



0 

1 

1 

1 



0 

2 

0 

2 



0 

2 

1 

0 



1 

0 

0 

2 


2 

0 

0 

2 

2 



0 

1 

1 

3 



0 

1 

2 

1 



0 

2 

0 

4 



0 

2 

1 

2 



0 

2 

2 

0 



1 

0 

1 

2 



0 

3 

0 

3 



0 

3 

1 

1 



1 

1 

0 

3 



1 

1 

1 

1 

3 

0 

0 

0 

0 

3 



0 

1 

0 

2 

1 

0 

0 

1 

CO 


0 

1 

0 

4 


16 



0 

0 

0 

1 

2 0 
0 
0 
0 
0 
0 
1 
0 
0 
0 
1 
1 


1 

2 

2 

0 

0 

1 

1 

2 

2 

2 

0 

3 

3 

3 

1 

1 


1 

0 

1 

0 

2 

1 

2 

0 

1 

2 

1 

0 

1 

2 

0 

1 


2 

3 
1 

4 
2 

5 
3 
1 

3 

4 
2 
0 
4 
2 
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TABLE II. The total binding energy and the distance of the e“e“e^ system. KL is the number 
of basis functions used. See text for the cases (i), (ii), and (iii) of basis functions. Atomic units are 
used. 


basis function 

K 

energy 


< r2 >1/2 



< r > 

radius 


e e 

e e"*" 

e~e~ 

(i) 

50 

0.26186974 







100 

0.26188326 

4.596 

6.968 

9.644 

5.494 

8.539 


150 

0.26188401 







200 

0.26188445 

4.597 

6.96882 

9.64479 

5.49475 

8.54000 

(ii) 

50 

0.26199779 







100 

0.26200455 

4.594 

6.958 

9.652 

5.489 

8.548 


150 

0.26200489 







200 

0.26200494 

4.595 

6.95835 

9.65285 

5.48963 

8.54856 

(iii) 

50 

0.26199953 







100 

0.26200465 

4.594 

6.958 

9.652 

5.489 

8.548 


150 

0.26200491 







200 

0.26200504 

4.595 

6.95812 

9.65254 

5.48957 

8.54846 


Hylleraas 

125 

0.262004895 

6.957 

9.650 

5.4891 

8.5476 

Hylleraas 

220 

0.2620050565 





CFHH II 

225 

0.262004673 

6.956 

9.650 

5.48881 

8.54699 

CE |26| 

800 

0.2620050702319 

6.95837 

9.65291 

5.4896332525 

8.5485806553 

Hylleraas |2|| 

744 

0.2620050702328 

6.95837 

9.65291 

5.489633252 

8.548580655 


18 













TABLE III. The total binding energies and the rms distances of the lowest S, P, D, and F 
states of the tt^ system. The mass set used is nit = 5496.918me, and = 206.7686me. Atomic 
units are used. 


L 

basis function 

JC 

energy 

radius 

< xlO^ 

tfi 

tt 

S 

(ii) 

200 

112.97253 

7.444 

11.15 

13.39 


(hi) 

200 

112.97300 

7.444 

11.15 

13.39 


CE in 

500 

112.9730179 




P 

(ii) 

200 

110.26189 

7.933 

11.82 

14.46 


(hi) 

200 

110.26210 

7.933 

11.82 

14.46 


CE m 

500 

110.2621165 




D 

(ii) 

200 

105.98288 

8.919 

13.18 

16.59 


(hi) 

200 

105.98301 

8.919 

13.18 

16.59 


(iii)“ 

200 

105.98292 

8.919 

13.18 

16.59 


CE 0“ 

2250 

105.982930 




F 

(ii) 

200 

101.43093 

10.71 

15.67 

20.39 


(hi) 

200 

101.43105 

10.71 

15.67 

20.38 


(ih)^ 

200 

101.43131 

10.71 

15.67 

20.39 


Adiabatic [0^ 


101.43 





“ mt = 5496.92158me,m^ = 206.768262me, 2i?oo = 27.2113961eV. 
^ rut = 5496.918me, = 206.769me, i?oo = 13.6058eV. 
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TABLE IV. The total binding energies and the rms distances of the lowest S, P, and D states of 
the td^ system. The mass set used is rtit = 5496.918me, = 3670.481me, and = 206.7686me. 

Atomic units are used. 


< 7-2 >1/2 xlO^ 


L 

basis function 

K 

energy 

radius 

dfj. 

tfj. 

td 

S 

(ii) 


100 

111.36357 

7.774 

11.73 

11.24 

13.92 





200 

111.36398 

7.774 

11.73 

11.24 

13.92 


(hi) 


100 

111.36363 

7.774 

11.73 

11.24 

13.92 





200 

111.36444 

7.774 

11.73 

11.24 

13.92 


GBCRC i 

1442 

111.364507 






CE 1 

13 


1400 

111.364511474 





P 

(ii 

) 


100 

108.17803 

8.417 

12.68 

12.03 

15.31 





200 

108.17914 

8.417 

12.68 

12.03 

15.31 


(hi) 


100 

108.17820 

8.416 

12.68 

12.03 

15.31 





200 

108.17940 

8.417 

12.68 

12.03 

15.31 


CE 1 

13 


1800 

108.1795424 






(iii 

)“ 


200 

108.17923 

8.417 

12.68 

12.03 

15.31 


GBCRC 11“ 

2662 

108.179385 





D 

(ii 

) 


100 

103.40632 

9.766 

14.81 

13.61 

18.19 





200 

103.40824 

9.769 

14.81 

13.61 

18.19 


(iii) 


100 

103.40733 

9.766 

14.81 

13.61 

18.19 





200 

103.40849 

9.769 

14.81 

13.61 

18.19 


GBCRC U 

1566 

103.408481 






mt = 5496.92158me,md = 3670.483014me, = 206.768262me, i?oo = 13.6056981eV. 
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TABLE V. The total binding energies and the rms distances of the lowest S, P, D, F, and G 
states of the helium atom (a e“e“). The correlated Gaussians of cases (ii) and (hi) are used for 
basis functions. The a-particle mass is = 7294.2618241me. Atomic units are used. 


< 7-2 > 1/2 


L 

basis function 

JC 

energy 

radius 

e a 

e e 

l^S 

(ii) 

200 

2.9033033 

0.8920 

1.093 

1.587 


(hi) 

200 

2.9033041 

0.8920 

1.093 

1.587 


CE 0 

300 

2.903304555 




2 ^S 

(ii) 

200 

2.1749299 

2.765 

3.386 

4.801 


(hi) 

200 

2.1749299 

2.765 

3.386 

4.801 


CE H 

over 350 

2.174930189 




2 ^P 

(ii) 

200 

2.1235432 

3.242 

3.971 

5.622 


(hi) 

200 

2.1235446 

3.242 

3.971 

5.622 


CE II 

over 350 

2.123545653 




2 ^P 

(ii) 

200 

2.1328785 

2.968 

3.635 

5.162 


(hi) 

200 

2.1328798 

2.967 

3.635 

5.162 


CE H 

over 350 

2.132880641 




3^D 

(ii) 

200 

2.0553385 

6.489 

7.949 

11.244 


(hi) 

200 

2.0553377 

6.489 

7.948 

11.243 


(ii)“ 

200 

2.0556201 

6.488 

7.948 

11.242 


(iii)“ 

200 

2.0556195 

6.488 

7.947 

11.242 


Hylleraas 

438 

2.05562073279 




3^D 

(ii) 

200 

2.0553538 

6.486 

7.945 

11.239 


(hi) 

200 

2.0553531 

6.486 

7.945 

11.238 


(ii)“ 

200 

2.0556355 

6.485 

7.944 

11.237 


(iii)“ 

200 

2.0556349 

6.485 

7.944 

11.237 


Hylleraas 

393 

2.05563630941 




4ii? 

(ii) 

200 

2.03097661 

10.963 

13.429 

18.992 


(hi) 

200 

2.03097596 

10.962 

13.427 

18.990 


(ii)“ 

200 

2.03125504 

10.962 

13.427 

18.989 


(iii)“ 

200 

2.03125314 

10.960 

13.425 

18.987 


Hylleraas 

438 

2.03125514434 




4^F 

(ii) 

200 

2.03097664 

10.963 

13.429 

18.992 
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(hi) 

200 

2.03097598 

10.962 

13.427 

18.990 


(h)“ 

200 

2.03125506 

10.961 

13.427 

18.989 


(iii)“ 

200 

2.03125317 

10.960 

13.425 

18.987 


Hylleraas |j^“ 

438 

2.03125516836 





(h) 

200 

2.0197237802 

16.589 

20.321 

28.738 


(hi) 

200 

2.0197114443 

16.553 

20.276 

28.675 


(h) 

200 

2.0197237803 

16.589 

20.321 

28.738 


(hi) 

200 

2.0197114445 

16.553 

20.276 

28.675 


“ The mass of the a-particle is assumed to be inhnite. 
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